Econometría II: ARIMA

Departamento de Economía

Carlos A. Yanes G.

2023-09-06

Paquetes con que se trabaja la sesión

R libraries

Los paquetes que vamos a utilizar y que se deben instalar para ser usados en clase son:

install.packages(c("pacman","TSstudio", "urca", "forecast", "devtools"))

# Si tiene problemas con TSstudio podemos mirar entonces en Github
devtools::install_github("RamiKrispin/TSstudio")

Ejecución

No olvide que debe desde luego tener presente en su documento de R Markdown cargar estos paquetes.

library(pacman)
p_load(tidyverse, TSstudio, tseries, urca, forecast)

Preambulo

Modelar ARIMA

  • Debemos primero mirar la estacionariedad.
  • Si no es estacionaria la serie hay que diferenciar.
  • Despues de testear procedemos a leer sus funciones de autocovarianzas denominadas correlogramas
  • Constatar el RUIDO BLANCO
  • Predecir

Datos

Cartera comercial de bancos

Agradecimientos

Esta sesión se construye a partir de los elementos de los paquetes de tseries y (Krispin 2019), quien desarrolló el paquete TSstudio respectivamente

A lo que vinimos

Cargar datos

library(pacman)
p_load(readxl, TSstudio, tidyverse, stats, urca, forecast, ggfortify, ggplot2, tseries)
bd <- read_excel("cartera.xlsx")
bd <- bd |> select(Cartera)
cartera <- ts(bd, frequency=12, start=c(2016,1))

Graficar normal

Code
ts.plot(cartera,
        title = "Cartera comercial de Bancos", 
        Ytitle = "En Miles de Millones de $",
        Xtitle = "Años")

Prueba de estacionariedad (1)

Code
prim <- ur.df(y=cartera, lags=0, type='none')
summary(prim)

############################################### 
# Augmented Dickey-Fuller Test Unit Root Test # 
############################################### 

Test regression none 


Call:
lm(formula = z.diff ~ z.lag.1 - 1)

Residuals:
   Min     1Q Median     3Q    Max 
 -8614  -1856   -262   1731  14834 

Coefficients:
         Estimate Std. Error t value Pr(>|t|)    
z.lag.1 0.0068402  0.0008212   8.329 1.14e-12 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 3653 on 86 degrees of freedom
Multiple R-squared:  0.4465,    Adjusted R-squared:  0.4401 
F-statistic: 69.38 on 1 and 86 DF,  p-value: 1.139e-12


Value of test-statistic is: 8.3292 

Critical values for test statistics: 
     1pct  5pct 10pct
tau1 -2.6 -1.95 -1.61
Code
prim2 <- ur.df(y=cartera, lags=0, type='drift')
summary(prim2)

############################################### 
# Augmented Dickey-Fuller Test Unit Root Test # 
############################################### 

Test regression drift 


Call:
lm(formula = z.diff ~ z.lag.1 + 1)

Residuals:
    Min      1Q  Median      3Q     Max 
-8991.9 -1598.5   365.7  1936.5 14876.5 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)   
(Intercept) -4.283e+03  2.453e+03  -1.746  0.08443 . 
z.lag.1      1.571e-02  5.144e-03   3.054  0.00302 **
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 3610 on 85 degrees of freedom
Multiple R-squared:  0.09887,   Adjusted R-squared:  0.08827 
F-statistic: 9.326 on 1 and 85 DF,  p-value: 0.003017


Value of test-statistic is: 3.0538 37.0386 

Critical values for test statistics: 
      1pct  5pct 10pct
tau2 -3.51 -2.89 -2.58
phi1  6.70  4.71  3.86
Code
prim3 <- ur.df(y=cartera, lags=0, type='trend')
summary(prim3)

############################################### 
# Augmented Dickey-Fuller Test Unit Root Test # 
############################################### 

Test regression trend 


Call:
lm(formula = z.diff ~ z.lag.1 + 1 + tt)

Residuals:
    Min      1Q  Median      3Q     Max 
-8860.2 -1405.5   318.1  1848.6 14698.6 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)
(Intercept) 8.965e+02  6.550e+03   0.137    0.891
z.lag.1     1.964e-04  1.890e-02   0.010    0.992
tt          4.831e+01  5.663e+01   0.853    0.396

Residual standard error: 3616 on 84 degrees of freedom
Multiple R-squared:  0.1066,    Adjusted R-squared:  0.08534 
F-statistic: 5.012 on 2 and 84 DF,  p-value: 0.008785


Value of test-statistic is: 0.0104 24.856 5.0119 

Critical values for test statistics: 
      1pct  5pct 10pct
tau3 -4.04 -3.45 -3.15
phi2  6.50  4.88  4.16
phi3  8.73  6.49  5.47
Code
# adf.test(cartera,k=0)
adf.test(cartera) # automático

    Augmented Dickey-Fuller Test

data:  cartera
Dickey-Fuller = -1.6971, Lag order = 4, p-value = 0.7006
alternative hypothesis: stationary

Analisis de funciones de autocorrelación

  • Permiten mirar o -sospechar- si la serie es estacionaria
  • La serie de cartera por ninguna parte parece ser estacionaria
Code
acf(cartera, lag.max = 20)

Code
pacf(cartera, lag.max = 20)

Comportamiento Lags

Code
ts_lags(cartera)

Comportamiento Lags

Code
ts_lags(cartera, lags=c(1,4,6,11))

Diferenciemos

Code
ts_plot(diff(cartera, lag = 1),
        title = "Transformación en Nivel",
        Xtitle = "Años",
        Ytitle = "Diferencia de cartera")

Diferenciación opcional

Code
ts_plot(diff(log(cartera), lag = 1),
        title = "Transformación en logaritmo",
        Xtitle = "Años",
        Ytitle = "Diferencia/logaritmica de cartera")

Nuevamente test e identificación

Code
dcartera<-diff(cartera, lag=1)
adf.test(dcartera)

    Augmented Dickey-Fuller Test

data:  dcartera
Dickey-Fuller = -3.0574, Lag order = 4, p-value = 0.1415
alternative hypothesis: stationary

############################################### 
# Augmented Dickey-Fuller Test Unit Root Test # 
############################################### 

Test regression trend 


Call:
lm(formula = z.diff ~ z.lag.1 + 1 + tt)

Residuals:
    Min      1Q  Median      3Q     Max 
-8738.8 -1752.5   150.9  1714.3 13670.4 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 512.7969   731.5826   0.701   0.4853    
z.lag.1      -0.6126     0.1010  -6.066 3.74e-08 ***
tt           31.7379    15.3805   2.064   0.0422 *  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 3336 on 83 degrees of freedom
Multiple R-squared:  0.3072,    Adjusted R-squared:  0.2905 
F-statistic:  18.4 on 2 and 83 DF,  p-value: 2.433e-07


Value of test-statistic is: -6.0661 12.2669 18.3991 

Critical values for test statistics: 
      1pct  5pct 10pct
tau3 -4.04 -3.45 -3.15
phi2  6.50  4.88  4.16
phi3  8.73  6.49  5.47

Identificación para pronostico

Code
par(mfrow=c(1,2)) 
acf(dcartera)
pacf(dcartera)

Modelos

Code
arima_1 <- arima(cartera, order = c(1,1,0))
arima_1

Call:
arima(x = cartera, order = c(1, 1, 0))

Coefficients:
         ar1
      0.6728
s.e.  0.0782

sigma^2 estimated as 12927179:  log likelihood = -836.05,  aic = 1676.11
Code
arima_2 <- arima(cartera, order = c(2,1,0))
arima_2

Call:
arima(x = cartera, order = c(2, 1, 0))

Coefficients:
         ar1     ar2
      0.4983  0.2529
s.e.  0.1041  0.1040

sigma^2 estimated as 12090598:  log likelihood = -833.2,  aic = 1672.41

Tenemos modelo!!

Test de Ruido Blanco

Code
forecast::checkresiduals(arima(cartera, order = c(1, 1, 0), include.mean = FALSE))

    Ljung-Box test

data:  Residuals from ARIMA(1,1,0)
Q* = 37.168, df = 17, p-value = 0.003193

Model df: 1.   Total lags used: 18

Pronostico de cartera

Code
arima_fc <- forecast(arima_1, h = 12)
plot_forecast(arima_fc,
 title = "Pronostico modelo AR(1)", 
 Ytitle = "Valores",
 Xtitle = "Años")

Detras de la serie

Code
# Crear el gráfico de las dos series de tiempo
plot.ts(cartera, col = "blue", xlab = "Tiempo", ylab = "Cartera Comercial", main = "Predicción de cartera comercial")
lines(arima_fc$fitted, col = "red")
legend("topleft", legend = c("Serie Original", "Predicción"), col = c("blue", "red"), lty = 1)

Valores

Code
arima_fc
         Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
May 2023       640573.6 635965.9 645181.4 633526.7 647620.6
Jun 2023       641772.5 632792.5 650752.5 628038.7 655506.2
Jul 2023       642579.0 629291.8 655866.3 622258.0 662900.1
Aug 2023       643121.7 625746.0 660497.3 616547.9 669695.4
Sep 2023       643486.8 622289.9 664683.6 611069.0 675904.5
Oct 2023       643732.4 618983.1 668481.6 605881.7 681583.1
Nov 2023       643897.6 615846.9 671948.4 600997.7 686797.5
Dec 2023       644008.8 612882.6 675134.9 596405.5 691612.1
Jan 2024       644083.6 610082.3 678084.8 592083.1 696084.0
Feb 2024       644133.9 607433.6 680834.2 588005.6 700262.2
Mar 2024       644167.8 604922.8 683412.7 584147.8 704187.7
Apr 2024       644190.5 602536.5 685844.6 580486.2 707894.9

Selección de modelos

  • Existen varias formas de buscar elegir un modelo entre varios que se han estimado por pensar que es nuestro mejor modelo generador de datos. Estos vienen a ser los criterios de Akaike (AIC) y Bayes (BIC)
  • \(AIC=-2\times log(L)+ 2\times p\)
  • P es el número de parámetros del modelo
  • L la razón de máxima verosimilitud del modelo
  • \(BIC=-2\times log(L)+ p \times log(n)\)
  • n es es tamaño de la muestra

Gracias por su atención!!

Referencias

Krispin, Rami. 2019. Hands-on Time Series Analysis with r: Perform Time Series Analysis and Forecasting Using r. Packt Publishing Ltd.